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Abstract: A simple method is proposed for determining the masses of new particles in 
collider events containing a pair of decay chains (not necessarily identical) of the form 
Z->Y + l,Y->X + 2, X->N + 3, where 1,2 and 3 are visible but N is not. Initial 
study of a possible supersymmetric case suggests that the method can determine the four 
unknown masses in effectively identical chains with good accuracy from samples of a few 
tens of events. 
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1. Introduction 

The discovery of new physics beyond the Standard Model is a primary objective of ex- 
periments at the Large Hadron Collider and other future colliders. In many models of 
BSM physics, a rich spectrum of new particles is predicted in the mass range accessible 
at the LHC. Many of these particles are weakly interacting and have quantum numbers 
that favour sequential decays into chains of other new particles plus visible jets and/or 
leptons. Typically the endpoint of the chain is a stable invisible particle that is a dark 
matter candidate. Classic examples are the squark decay chain in supersymmetric models, 

q^xl + Q, xl^P + F, £ t -X? + * ± , (1-1) 

where the neutralino Xi is the lightest supersymmetric particle (LSP), and the excited 
quark decay in models with universal extra dimensions, 

q* -> Z* + q , Z* — ► (*~ n + , r± 7* + £*= , (1.2) 

where the photon excitation 7* is the lightest Kaluza-Klein particle (LKP). 

Determining the masses of the new particles in such decay chains, especially the dark 
matter candidate, is clearly of great importance. Many approaches to this problem have 
been proposed [1-47], based mainly on the measurement of endpoints or other features in 
the distributions of invariant masses or specially constructed observables, or on explicit 
solution for the unknown masses using multiple events. 

The present paper investigates a somewhat different approach which is particularly 
suited to processes in which there are two three-step decay chains of the form ( |1.1| ) or ( |1.2j ) . 
In principle the chains need not be the same, although in practice it would be too much to 
expect to determine the eight masses involved in non-identical chains. We shall see that 
for practically identical chains, such as first- and second-generation squark pair production 
and decay as in (|1 . 1| ) , determination of the four sparticle masses and reconstruction of the 
LSP momenta appears possible with reasonable numbers of events. 
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Figure 1: Event topology 



2. Method 

Consider the double decay chain in fig. H. The 4-momenta in the upper chain should satisfy 

(Pi + P2 + P3 + Pif = M\ 
{P2+PZ+PA? = My 

{pz+Pi? = M\ 



p\ = Ml 



(2.1) 



Leaving aside the last equation, the others give three linear constraints on the invisible 
4-momentum p±: 

- 2p x • p 4 = M Y -M 2 Z + 2pi • p 2 + 2pi • p 3 + m? = S x 
-2p 2 -Pi = M 2 X -M Y + 2p 2 ■p 3 + ml = S 2 



-2p 3 ■ p 4 = M 2 N -M 2 x +ml = S z 



(2.2) 



Similarly for the lower chain 

- 2p 5 ■ p 8 = My, - Affy + 2p 5 ■ p 6 + 2p 5 ■ p 7 + ml = S 5 
-2p 6 ■ p 8 = M\, - M Y , + 2p 6 ■ p 7 + wig = 5 6 
-2p 7 ■ p 8 = M%, - M\, + m 2 = S 7 

We also have the missing transverse momentum constraints 



(2.3) 



P4 + P8 — Pmiss 

y , y y 

Pa ' P% — Pmiss 



Si 



Let us make an 8-vector of the invisible 4-momenta, 

P = (p X 4,Pa,P Z 4, Ea,p^,pI,p z 8 , E 8 



(2.4) 



(2.5) 
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Then we have 



where A is the 8x8 matrix 



A = 2 
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p = A _1 S = DM + E 



(2.6) 
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Furthermore S may be written as 

S = BM + C 

where M is the vector of masses-squared to be determined, 

M = (M 2 ,M 2 ,M 2 X ,M 2 N ,M 2 ,,M 2 ,,M 2 ( ,,M 2 N ,) , 
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and 

C = (2pi ■ p 2 + 2pi ■ p 3 + ml, 2p 2 ■ P3 + m 2 ,, ml,p^ iss , 
2p 5 ■ p 6 + 2p 5 • p 7 + mj, 2p 6 ■ p 7 + ml, mf,^ iss ) . 

Hence the solution for the invisible 4-momenta is 



where D = A^B and E = A _1 C. 

Since A n for each event n is a sparse matrix it is easy to invert it and hence to obtain 
and store D„ and E n (72 real numbers) for each event. Then for every hypothesis for the 
unknown masses M we immediately get a unique solution for the invisible 4-momenta in 
each event, P n . For the correct mass hypothesis, these satisfy the quadratic constraints 

(Pl)n = (Pi ~ Pi ~ Pi - Pi)n = M 2 N 

(Pl)n = (Pi ~ Pi - Pi - Pl)n = M 2 N , (2.13) 
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We can therefore measure the goodness of fit for the mass hypothesis M by the quantity 1 
e(M) = £ [(pg) B - M%] 2 + £ [(pg) B - M 2 N> ] 2 (2.14) 

n n 

The method is then to find the best-fit hypothesis for the masses by minimizing £ 2 . Note 
that this quantity tests the goodness of fit to all the masses equally, since for example it 
follows from eq. Q2.2j) that 



(P3 + Pi)l ~M 2 X = (pl)n - M 2 N (2.15) 

and similarly for all the other unknown masses. 

To account for combinatorial ambiguities, we must evaluate £ 2 for all permitted particle 
combinations for each event. Notice however that a different combination within one chain 
corresponds to a permutation of the rows of the matrix A. Therefore the inverse is given 
by the same permutation of the columns of A" 1 , and no extra inversions or matrix storage 
are necessary. In the case that the mass difference between particles 2 and 3 is negligible 
(as for dileptons), the vector C is invariant under their exchange; similarly for 6 and 7. 
Thus only combinations that exchange particles between chains require extra data storage. 

For every mass hypothesis, we should use the lowest value of £ 2 amongst all the allowed 
combinations for each event. At the best-fit point, this also shows which combination is 
most likely to be the correct one. Then the corresponding reconstructed momenta can be 
used, for example, to test spin hypotheses for the particles involved. 

This method is closest in spirit to that of refs. [28,44], in which pairs of events with 
the same decay chains are used to solve explicitly for the unknown masses. As this is a 
non-linear problem, each pair yields multiple solutions, which are narrowed down to the 
correct one as more pairs are solved. The advantage of the present method is that each 



event contributes independently and additively to the goodness-of-fit function (2.14), which 
is obtained by a simple linear computation for any mass hypothesis. The problem then 
reduces to the familiar one of function minimization. The additivity property also means 
that combining event samples and the statistical interpretation of results become more 
straightforward. Furthermore the results from other methods such as edge analyses can 
easily be included as constraints on the mass parameters in the minimization search. 

3. Results 

As an illustration of the method, it was applied to the process of squark-pair production 
at the LHC (pp at 14 TeV centre-of-mass energy). The SUSY mass spectrum and decay 
branching ratios were taken to be those of CMSSM point SPS la [48]. At this point the 
SUSY production cross section at the LHC is about 50 pb and there is a good probability for 
squark production and decay into quark jets and dileptons via ( |1.1| ). Some of the squarks 



1 The symbol f 2 rather than x 2 reminds us that it has no probabilistic interpretation. However, in 
an experimental analysis event-to-event variations in momentum uncertainties could be propagated from 



eqs. (2.7) and (2.11) into (2.14) to give greater weight to events with higher precision, thus defining a 



quantity more like x ■ My thanks to Ben Gripaios for this suggestion. 
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Figure 2: Best-fit results for SUSY point SPS la. Each point represents a sample of 25 events. 
The dashed lines show the true masses. 



are produced directly and some come from gluino decay; the production mechanism affects 
their momentum and rapidity distributions but is otherwise irrelevant for our purposes. 
Decays of the two squarks into unlike dileptons (e + e~ fi + fi~) were selected to limit the 
number of allowed combinations of jets and near and far leptons to eight, necessitating the 
storage of 144 real numbers for each event as explained above. 

Third-generation squarks were excluded, as their different masses prevent a good fit 
with a single squark mass. Experimentally, this would involve vetoing events with a tagged 
6-jet. At SUSY point SPS la only left-squarks have significant branching ratios into the 
mode (1.1) and so the left-right squark mass splitting is not a problem here. Therefore a 
four-parameter fit with Mj + 4 = Mi in eq. ( |2.9| ) is appropriate. 

Figure [2] shows the best-fit results for 100 Monte Carlo samples of 25 events each, 
selected as described above. A sample of 25 such events corresponds to an integrated 
luminosity of about 3 fb _1 . The events were generated with HERWIG version 6.510 [49-51] 
and the jet and lepton momenta used in the analysis were at parton level (after parton 
showering but before hadronization), with perfect jet reconstruction and no momentum 
smearing in this figure. The missing transverse momentum was taken to be that of the 
LSPs alone, again without smearing. However, HERWIG smears all unstable particle masses 
with the appropriate Breit-Wigner distributions, an effect that is significant for the squarks 
(Tg ~ 5 GeV) but negligible for the sleptons and neutralinos. 

We see that in this idealized situation the four new particle masses are usually quite well 
determined. As summarized in the first row of table |l], the r.m.s. variation in the estimated 
mass is 20 GeV for the squark and around 10 GeV for the slepton and neutralinos, with 
mean values within 1 or 2 GeV of the true ones. The best-fit combination is the correct 
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Table 1: Fitted masses and r.m.s. variations for samples of 25 events. The true average masses 
are shown in the heading (all in GeV) . The quantity is the fraction of samples surviving the £ 2 
cut, while / cor is the fraction of events with the correct best-fit combination. 

one (at the best-fit point) in 72% of events. The fits are not perfect, and there are incorrect 
choices of combination, because of the intrinsic differences amongst the squark and slepton 
masses and the Breit-Wigner smearing. 

The search for the best fit is rather tricky because the ^-surface is not smooth, owing 
to sudden changes in the best-fit combinations as the mass parameters are varied. For this 
study, the SIMPLEX method in MINUIT [52] was used. The surface would be smooth if 
one added the £ 2 contributions of all combinations, but then the sensitivity to the correct 
solution is reduced and biases are introduced by the huge contributions of wrong combina- 
tions. 

The mass resolution can be improved by eliminating data sets with large best-fit values 
of £ 2 . For example, if we require our 25 event sample to have total £ 2 < 100 in units of 
(100 GeV) 4 , then 80% of the samples survive and the fluctuations in the fitted masses are 
reduced as indicated in the second row of table [l]. 

As a rough indication of the possible effects of hadronization, reconstruction errors 
and detector resolution, the jet and lepton momenta and also the missing momenta were 
smeared with a gaussian distribution of r.m.s. width 5p/p = 5% or 10%, with the results 
shown in the lower rows of table [l] and, for 10% smearing, in fig. |3[ 

The fluctuations in the fitted masses naturally increase with increasing smearing, and 
the fraction of correct best-fit combinations decreases. There is a slight bias of the mass 
estimates towards lower values due to momentum smearing. As in the unsmeared case, 
a cut on £ 2 reduces the fluctuations at the expense of rejecting a fraction of the event 
samples. This means, for example, that in a single experiment with 25 events and 10% 
smearing there is about a 58% chance that £ 2 > 200, in which case the estimates of the 
masses are likely to be poorer than those shown in the last row of the table. 

Since the individual events are statistically independent, the fluctuations in the best-fit 
mass estimates decrease inversely as the square root of the size of the event sample. Any 
systematic biases in these estimates would not decrease with improving statistics, but the 
results with momentum smearing suggest that such effects should be small. They could 
be corrected using more detailed Monte Carlo simulations if present in a real experimental 
analysis. 
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Figure 3: As in figure 2, but with momentum smearing Sp/p = 10%. 

One possible source of bias would be the neglect of jet invariant mass in the recon- 
struction of the quark jets. In the present study, the effect of this was investigated by 
rescaling either the jet momentum or the jet energy to give zero jet mass. In the former 
case there was little effect, but energy rescaling resulted in a downward bias of about 25 
GeV in the estimated squark mass and around 5 GeV in the other masses. 



4. Conclusions 

The method of mass determination presented above is simple to apply and looks promising 
for the class of processes studied here. It can readily be extended to more complicated 
final states involving different or longer decay chains. Combinatorial background does not 
appear to be a serious problem although other backgrounds and the effects of additional 
jets due to QCD radiation remain to be investigated. These would be best studied in the 
framework of full simulations including detector effects. Comparative studies along these 
lines for this and other mass determination methods are in progress [53]. 

A variant of this method can also be applied to events with one three-step decay chain 
like ( [1.1] ) or (1.2) and one shorter (two-step) chain, for example 



9^q' + q', q'^Xi + q', (4.1) 

or analogously 



+ q , q -> 7 +q . (4.2) 

In this case the four constraints (|2.1| ) on the longer chain can be solved for the invisible 
4-momentum p4, with a two- fold ambiguity since one constraint is now quadratic. The 
kinematics of the shorter chain can then be reconstructed and the goodness of fit ( 2.1 -j| ) 
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computed for both solutions. Choosing the solution with the better fit for each mass 
hypothesis, one can now proceed as in the case of two three-step chains. Further discussion 
of this case will be presented in a later paper [54]. 
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